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Fusion of Sparse Reconstruction Algorithms 
for Multiple Measurement Vectors 

Deepa K G, Sooraj K. Ambat, and K.V.S. Hari 


Abstract 

We consider the recovery of sparse signals that share a common support from multiple measurement 
vectors. The performance of several algorithms developed for this task depends on parameters like 
dimension of the sparse signal, dimension of measurement vector, sparsity level, measurement noise. 

We propose a fusion framework, where several multiple measurement vector reconstruction algorithms 
participate and the final signal estimate is obtained by combining the signal estimates of the participat¬ 
ing algorithms. We present the conditions for achieving a better reconstruction performance than the 
participating algorithms. Numerical simulations demonstrate that the proposed fusion algorithm often 
performs better than the participating algorithms. 

Index Terms 

Compressed Sensing, Fusion, Sparse Signal Reconstruction, Multiple Measurement Vectors 

I. Introduction 

Consider the standard Compressed Sensing (CS) measurement setup where a A-sparse signal 
x G M Vx 1 is acquired through M linear measurements via 

b = Ax + w, (1) 

where A G R Mx ,v denotes the measurement matrix, b G M Mx 1 represents the measurement 
vector, and w G M Mxl denotes the additive measurement noise present in the system. The 
reconstruction problem, estimating x from ([B using A and b, is known as Single Measurement 
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Vector (SMV) problem. In this work, we consider the Multiple Measurement Vector (MMV) 
problem flU where we have L measurements: b^ 1 ' = Ax^ + w (1 \ t>( 2) = Ax (2 ^ + w (2 \ • • •, 
b^ = Ax (L ) + w {L \ The vectors {x ( ^}^ 1 are assumed to have a common sparse support- 
set. The problem is to estimate x® (7 = 1,2Instead of recovering the L signals 
individually, the attempt in the MMV problem is to simultaneously recover all the L signals. 
MMV problem arises in many applications such as the neuromagnetic inverse problem in 
Magnetoencephalography (a modality for imaging the brain) Q, Q, array processing H, 
non-parametric spectrum analysis of time series O, and equalization of sparse communication 
channels m. 

Recently many algorithms have been proposed to recover signal vectors with a common sparse 
support. Some among them are algorithms based on diversity minimization methods like /: 2 . i 
minimization If7l . and M-FOCUSS 0}, greedy methods like M-OMP and M-ORMP [U, and 
Bayesian methods like MSBL (§1 and T-MSBL [0. 

However it has been observed that the performance of many algorithms depends on many 
parameters like the dimension of the measurement vector, the sparsity level, the statistical 
distribution of the non-zero elements of the signal, the measurement noise power etc. fl9]|. Thus it 
becomes difficult to choose the best sparse reconstruction algorithm without a priori knowledge 
about these parameters. 

Suppose we have the sparse signal estimates given by various algorithms. It may be possible 
to merge these estimates to form a more accurate estimate of the original. This idea of fusion of 
multiple estimators has been proposed in the context of signal denoising in lITOll where fusion was 
performed by plain averaging. Recently, Ambat et al. lUTll — [TT31l proposed fusion of the estimates 
of sparse reconstruction algorithms to improve the sparse signal reconstruction performance of 
SMV problem. 

In this paper, we propose a framework which uses several MMV reconstruction algorithms and 
combines their sparse signal support estimates to determine the final signal estimate. We refer to 
this scheme as MMV-Fusion of Algorithms for Compressed Sensing (MMV-FACS). We present 
an upper bound on the reconstruction error by MMV-FACS. We also present a sufficient condition 
for achieving a better reconstruction performance than any participating algorithm. By Monte- 
Carlo simulations we show that fusion of viable algorithms leads to improved reconstruction 
performance for the MMV problem. 
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Notations: 

Matrices and vectors are denoted by bold upper case and bold lower case letters respectively. 
Sets are represented by upper case Greek alphabets and calligraphic letters. A r denotes the 
column sub-matrix of A where the indices of the columns are the elements of the set T. X 7 -. 
denotes the sub-matrix formed by those rows of X whose indices are listed in the set T. X A 
is the matrix obtained from X by keeping its K rows with the largest f 2 - n orm and by setting 
all other rows to zero, breaking ties lexicographically. supp(X) denotes the set of indices of 
non-zero rows of X. For a matrix X, x !/l denotes the I th column vector of X. X, denotes the 
reconstructed matrix by the z lh participating algorithm. The complement of the set T with respect 
to the set { 1 , 2 ,..., N} is denoted by T c . For two sets T\ and T2, T\ \ T2 = 7{ D T£ denotes 
the set difference. \T\ denotes the cardinality of set T. A^ and A T denote the pseudo-inverse 
and transpose of matrix A, respectively. The ( p , q) mixed norm of the matrix X is defined as 



The Frobenius norm of matrix A is denoted as ||A|| F . 

II. PROBLEM FORMULATION 

The MMV problem involves solving the following L under-determined systems of linear 
equations 

b (/) = Ax w + w (/) , 1 = 1,2,3,..., L (2) 

where A e R MxiV (M <R N ) represents the measurement matrix, b- G R Mx 1 represents the I th 
measurement vector, and x ll> G R JVxl denotes the corresponding A-sparsc source vector. That 
is, |supp(x^)| < K and share a common support-set for l = 1,2w ( ^ G R Mxl 
represents the additive measurement noise. We can rewrite (Eb as 

B = AX + W (3) 

where X= [x^, x( 2 ),..., x^], W = [w^, w< 2 ),..., w^], and B = [b«, b( 2 ),..., b^]. For 
a matrix X, we define 

L 

supp(X) = |^J supp(x*). 

i— 1 

In ©, we assume that X is jointly K -sparse. That is, |supp(X)| < K. There are at most K 
rows in X that contain non-zero elements. We assume that K < M and K is known a priori. 
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III. Fusion of Algorithms for Multiple Measurement Vector Problem 

In this paper, we propose to employ multiple sparse reconstructions algorithms independently 
for estimating X from © and fuse the resultant estimates to yield a better sparse signal estimate. 
Let P > 2 denote the number of different participating algorithms employed to estimate the 
sparse signal. Let % denote the support-set estimated by the z lh participating algorithm and let T 
denote the true-support-set. Denote the union of the estimated support-sets as T, i.e., T = U[ =1 71, 
assume that R = F < M. We hope that different participating algorithms work on different 
principles and the support-set estimated by each participating algorithm includes a partially 
correct information about the true support-set T. It may be also observed that the union of the 
estimated support-sets, T, is richer in terms of the true atoms as compared to the support-set 
estimated by any participating algorithm. Also note that, once the support-set is estimated, the 
non-zero magnitudes of X can be estimated by solving a Least-Squares (LS) problem on an 
over-determined system of linear equations. Hence if we can identify all the true atoms included 
in the joint support-set T, we can achieve a better sparse signal estimate. 

Since we are estimating the support atoms only from F, we need to only solve the following 
problem which is lower dimensional as compared to the original problem ©: 

B = ArXr, + W, (4) 

where Ap denotes the sub-matrix formed by the columns of A whose indices are listed in 
T, Xp . denotes the submatrix formed by the rows of X whose indices are listed in F, and 
W = W + ApcXpc.. The matrix equation © represents a system of L linear equations which 
are over-determined in nature. We use the method of LS to find an approximate solution to the 
overdetermined system of equations in ©. Let Vp : denote the LS solution of ©. We choose 
the support-set estimate of MMV-FACS as the support of V A , i.e., indices of those rows having 
the largest A-norm. Once the non-zero rows are identified, solving the resultant overdetermined 
solution using LS we can estimate the non-zero entries of X. MMV-FACS is summarized in 
Algorithm [0 
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Algorithm 1 : MMV-FACS _ 

Inputs: A e R MxN , B e M MxL , I<, and 

l ) i=l:P 

Assumption: | U %\ < M. 

i =1 

Initialization: V = 0 e R Nxl . 

Fusion: 

1: r = u t; 

2—1 

2: Vr,; = Aj.B, V r c, ; = 0; 

3: T = supp(V A ); 

Outputs: T and X (where X r . = At. B and X, p . . = 0) 


Remark: 


An alternate approach for solving an MMV problem is to stack all the columns of B to get 
a single measurement vector. Then © in a noiseless case becomes 


bi 

b 2 


1 

> 

> 

1 

O 


Xl 

x 2 

_b l_ 

MLx 1 

0 

A 

A 

MLxNL 

X/._ 


NLxl 


where bj and x,- (i = 1,2 ,,L) denote the i lh column of B and X respectively. Now, we have 
the following SMV problem. 




" A 



' h' 


f\ 


Xi 



A U 



b 2 

~ 



X 2 



0 A 



pNL_ 




XNL_ 


MLx 1 

A 






MLxNL 



(5) 


J NLxl 


In principle, we can solve © using FACS with sparsity level LK. Note that, after stacking X 
column-wise, we lost the joint sparsity constraint imposed on X in the MMV problem in ©. The 
LK non-zero elements estimated from © using FACS can be from more than K different rows 
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of X. In the worst case, the estimate of FACS may include non-zero elements from min (LK, M ) 
different rows of X. Then we will end up with an estimate of X with LK non-zero rows, which 
is highly undesirable. Hence stacking the columns of the observation matrix B and solving it 
using FACS is not advisable. Note that Step 3 in Algorithm Q] ensures that MMV-FACS estimates 
only K non-zero rows of X. 


IV. Theoretical Studies of MMV-FACS 


In this section, we will theoretically analyse the performance of MMV-FACS. We consider 
the general case for an arbitrary signal matrix. We also study the average case performance of 
MMV-FACS subsequently. 

The performance analysis is characterized by SRER extended for MMV which is defined as 


SRER = 


IXI 


x-x 


2 1 


( 6 ) 


where X and X denote the actual and reconstructed signal matrix respectively. 


Lemma 1. Suppose that A satisfies the relation, for some constant 8r+k £ (0,1), 

||AX|| f < \A + 5r + k ||X || f , 


where ||X || 0 < R + K and 5r+k £ (0,1). Here ||X || 0 denotes the number of non-zero rows of 
the matrix X. Then, for every matrix X, 


AX|| f < y/l + 8r+k 



1 

VR+K 



2,1 


Proof : Proof is given in Appendix lAl 


Lemma 2. Consider A e M MxAr and let 71 & be two subsets of {1,2,.. .N} such that 
7 iD72 = 0. Assume that 5|Ti|+|75| < 1 and let Y be any matrix, such that span(Y ) G sparf A-j] ) 
and R = Y — A^A^Y. Then we have 


1 - 


V1I+IT2I \ 


1 “ Vd+ \T2\J 

Proof : Proof is given in Appendix lAl 


|Y|| 2 < 


|R|| 2 < l|Y|| 2 
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A. Performance Analysis for Arbitrary Signals under Measurement Perturbations 

We analyse the performance of MMV-FACS for arbitrary signals and give an upper bound on 
the reconstruction error in Theorem Q] We also derive a sufficient condition to get an improved 
performance of MMV-FACS scheme over any given participating algorithm. 


Theorem 1. Let X be an arbitrary signal with T = supp(X A ). Consider the MMV-FACS 
setup discussed in Section [777] and assume that the measurement matrix A satisfies RIP with RIC 
Sr+k- We have the following results: 

i) Upper bound on reconstruction error: We have, 


X-X 


< c x ||X - X A \\ F + C 2 ||X - X A || 2 x + C 3 ||Xrc, : || F + !/ II W|| 


where C x = (l + vyj\ + 6 R+K ^, C 2 = ? C: 


1 + 5 


R+K 


3-5 


(1 — 5 R+K ) 2 


, and v = 


R+K 


(1 — Sr+k) 2 

ii) SRER gain: 


For 


X 




f 0 and ||X r c : |L f 0, MMV-FACS provides at least SRER gain of 


f 1 S \ 2 \ 2 

7 ---- ] over the i th participating algorithm if 

(1 + 8r+k + 3C + 3pr]iJ 

Vi < 77-i- R+ n} — F 7 \’ where Vi — ~, C = 7777——, and 


(1 + $r+k + 3C + 3£) 


£ — (3a/1 + 5 r+k + 1 


X-X 


K I 


||X rC ,:| 

F r 

x tc . 

' i 5 " 

> s 

F 


X 


r c ,:\\p 


3 ||X 


+ 


JR+K 


X-X 


K\ 


12,1 


r c ,: 1 + 


VR + K ||x rv | 


Proof: 
i) We have, 


Consider, 


X-X 


< X-X 


K I 


+ 


x"-x 


(7) 


X K -X 

< 

F 

(X'+,:-Xf, : 

+ 

F 


< 

- X t , ; 

+ 

F 


x^) tC: -x tc>: 


X 


K\ 


>T C , : 


(yX tv = 0) 


( 8 ) 
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Using the relations X r . — At B (from Algorithm CD) and A1 A, — I, we get 




:x A ') t - v.ii 


t T 


(X\ : -A^(AX + W) 


B = AX + W) 


- At (AX k + A(X - X A ') + W) 


X A ) t: - A\ A r (X K ) t : + A tc (X A ) tc : + A(X-X A ) + W 


A^A tc (X^) tCj; + At_A(X - X A ) + A^W 


K\ 


< 


(AfA^y 1 A»A fc (X K ) f ^ : 


+ 


At_A(X - X 


K\ 


+ 


At_W 


(9) 


Let denote the z lh column of matrix X and w f '' denote the i th column of matrix W, i = 
1, 2,... L. Now from Proposition 3.1 and Corollary 3.3 of |[T6ll we obtain the following relations. 


A^w« 


< 


w 




A t A ( ; 


r(0 


^) K ) 


< 


R+K 

A fx (i) - (x 


w )rillt 


( A t A t) A " A f (( x<i) ) K )f» 

Consider (fTOl) . we get 


< 


\/l — £r+/i 
5.R+A' 


1-5 


a+a: 


||((x«)*) 


r c 11 2 


( 10 ) 

(ID 

( 12 ) 


At_ w 0) 


< 


w 


2 1 — 5r+k 

Summing the above equation over i — 1, 2,... L, we obtain 


Vi = l,2,...L 


X|l A t 


(0 


1=1 


W 


a(_w 


a(_w 


< 


1-5 


R+K 


X 

i=l 


W 


< 


< 


1 — $R+K 
1 


y/ 1 3 ^ 


l|W||5. 

= I|W|| F 


R+K 


Similarly, summing the relations in (fill) and we obtain 


A*_A (X - X A ') 


< 


|a(x-x a ' 




(13) 


(14) 


R+K 


April 8, 2015 


DRAFT 


















































































(15) 


(A^A t ) 1 A^A tc (X A ') tc>: 

Substituting (fT3l).(fl4l) and (fl5l) in (HI), we get 

5r+k 


< 


JR+K 


1-5 


R+K 


X 


K\ 


ix 


K\ 


>T,: 


X 


r,: 


< 

F 1 — 3r+k 
5r+k 


;x k ) 


T+: 


< 


1-5 


R+K 


:x 


K\ 


>T C , : 


+ 


+ 


|a(x-x a ') 

\/l — $ R+K 
1 


)T C , : 


+ 


IWI 


y/tt. 


R+K 


1-5 


R+K 


(||A(X-X* 


|W|| F ) 


(16) 


Substituting (fT6l) in ([8]), we get 

1 


X K -X 


< 


(X A ) 


F 1 — 5 r+k 

Next, we will find an upper bound for 


+ 


r ’' f 1 — 5 r+k 
(X A > C: " 


(||A(X-X a ')|| f + ||W|| f ) (17) 


Define 7 a — T \ T. That is, 7 a is the set formed by the atoms in T which are discarded by 
Algorithm [D Since T Cf, we have T c = T c U 7 a and hence we obtain 


fX 


K\ 


>T C ,: 


< (X 


K\ 


> • 11 F 


+ 


fx 


K\ 




(18) 


We also have, 



< 

F 

( Vr .=)t Al : 

+ 

F 


< 

( V v)f A , : 

+ 1 
F 


(Vr, - (X K )r,:) tAi: 

+ ||Vr, : -(X A ')r, ; || F 


(19) 


Note that (Vr, : )j- . contains the 77-rows of Vr, : with highest row f 2 -norm. Therefore, using 
\T\ = \T\= K, we get 


(Vr, : ) 


7a,: 


< 


F 

K\ 


(Vr,:) r \ Tj . 

= ||V r \r,: — (X A )r\r,:| 
< ||(Vr,:-(X K )r, : )| 


(X*)r\r,: = 0) 


Substituting (l20l) in (IT9l) . we get 


:x 


K\ 


>Ta,: 


— 2 || (V r , : — (X A ) r , ; ) || F . 


( 20 ) 


( 21 ) 
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Now, consider 

||V r , : -(X x ) r , : | 


aIb - (x 


K\ 


Jr,: 


A{,(AX + W)-(X i 


r K\ 


- K\ 


A{, (AX a + A(X - X A ') + W) - (X ;r, : 
Ap^A r (X A ) r . + A r c(X Ar ) r c ) . + A(X — X 
Af,A r c(X K )rc, : + A+A(X - X K ) + A*W 


(X A ')r,: 

, t A _ _ 


ApA r = I) 


< 


ApA r c(X A ) r c >: 


A[A(X - X 


K\ 


+ 


ApW 


Using (fl3l) . (fl4l) and (IT5T) in (l 22 l) . we get 

Sr+K 


Vr,:-(X")r,:L< 


< 


1 — Sr+k 

&R+K 


x )r«=,:|| F + 


X A ) 


l|A(X-X^)|| ; 

\/l ~ &R+K 

l|A(X-X^)|| j 


+ 


liwih 


v/w 


R+K 


1 _ &R+K 

0 < 1 — 6r+k < 1 ) 

Using (|2T1) and (1231) in (fl 8 l) . we get 


+ ^l|W|| F 


1-5 


R+K 


1-5 


R+K 


(X A ) 


T c ,: 


F 1 — Or + k 


If+ i_5 


R+K 


(||A(X-X k )|| f + ||W|| p ) 


Let denote the i to column of matrix X A . The, we have, 


K 


Ax 


(0 


< (1 + $r+k) 


r (i) 


(;.• A satisfies RIP) 


1=1 


eK ~ r.a + $r+k) 


i=i 




K 


12 
If 


IIAX A ||p, < (1 + 5r+k) ||X 

Using Lemma [Hand (l25l) . we get 

||A(X — X^)^ < yjl + 5r + k Ilx-X 


K\ 


r K 


F+ Vr + k'^ X Xj /i|2 ’1 


( 22 ) 


(23) 


(24) 


(25) 


(26) 
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Substituting (1241) 


x A -x 


4J) in (1171) . we get 

< * + i « +K , 2 ||(XV,:|| F + n 3 [||A(X 

(1 — 5r+k) 2 " F (1 — Sr+k) 2 l 

^ ^urw U V^ + ^A+A |X — i 


||x-x A || F 


- X L + 

+ H|w|| f 


Viul 


x - x * 1 k, 


(using (l26l) ) 


< V\/\ + <7r+at | 

1 + ||v M 
+ (i-<w l|Xrvllf " 

. 3 — 5r+k 

where v = 71-a -Ti¬ 

ll — Or+k) 

Substituting (l27l) in ([7]) and using the definitions of Ci,C 2 , and C 3 , we get 

< c x ||x - x A || F + c 2 ||x - x A || 21 + C 3 ||X r c i: || F + Z/ ||W|| F . 


(27) 


x-x 


(28) 


ii) Using (1281) and 
X-X 


the definitions of 7 and //,, we get 

1 + $r+k + 3C + 3£ M n 
-----—-ll^r c , : || F 


F (1 — &r+kY 
1 + &R+K + 3£ + 3£ 
(1 — Sr+k ) 2 

< 1 + $R+K + 3£ + 3£ 


(1 — $r+k) z 


(X-XOt* . 

(X - X 


(X 


■VT C ,: 


0 ) 


V ^ n^-ix j 

Hence, we obtain the relation for SRER for MMV-FACS, in case of arbitrary signals, 

I I -*7- I I 2 

SRER | mmv-facs 


as 


IIXI 


X-X 


> 


IIXI 


x-x 


F _ x 
2 x 


( (i-w ) 2 V 

V (1 + $r+k + 3<T + 3 


Hence MMV- 


-FACS provides at least SRER gain of (--— $ R+I <) — — 

\{1 + 5 r+k + 3C + 3^)i7i 

rithm if r/ ?; < - ^ ^ R+K } -—. 

(1 + $r+k + 3C + 3£) 

Note that 7 - (±Z*B±kL < l. 

(1 + $r+k + 3£ + 3£) 


over z 1 * 1 algo- 


□ 


B. Exactly K-sparse Matrix 

Theorem Q] considered the case when X is an arbitrary matrix. If X is a /v'-sparse matrix 
then we have X = X A and £ = 0. Thus, it follows from Theorem Q] that, MMV-FACS 
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provides at least SRER gain of 


(1 — Sr+kY 


(1 + &R+K + 3C )Vi 


over z lh participating algorithm if Vi < 


(1 ~ 3r+k) 2 
(1 + $r+k + 3£) 

z lh Algorithm for a A'-sparse matrix is greater than that of an arbitrary matrix by a factor of 


Thus, the improvement in the SRER gain provided by MMV-FACS over the 


1 + 


3? 


(1 + 8r+k + 3£) 


The second part of Theorem [Q considers the case when 


X 


TV 


7 ^ 0 and ||X r c ; .|| F f 0. If 


X 


T c • 

' i. i • 


7 ^ 0 , then 71 ^ T. Also, ||X r =, : || F = 0 implies T C T. Suppose 


X 


T. c ■ 


= 0 , then 


the support-set is correctly estimated by /' lh algorithm and further performance improvement is 
not possible by MMV-FACS. Hence we consider the case where ||X r c ; .|| F = 0, and derive the 
condition for exact reconstruction by MMV-FACS in the following proposition. 


Proposition 1. Assume that ||Xrc, : || F = 0 and all other conditions in Theorem \l\ hold good. 
Then, in clean measurement case (W = 0), MMV-FACS estimates the support-set correctly and 
provides exact reconstruction. 


Proof'. We have 


X rc , : =0 ^TcT (29) 

B = A r X t: + W (30) 

From Algorithm [U we have V e M JVxL where V rc . = 0, and 

V r , : = A ] r B 

= Af (A r X t: + W) (using m) 

= Ap(A r X r , : + W) (using ([291)) 

= X r ,: + ApW. 

If W = 0 , then Vr, : = Xr, : and V = X (•/ T C T). Thus MMV-FACS estimates the support-set 
correctly from V. □ 

In practice, the original signal is not known and hence it is not possible to evaluate the 

performance w.r.t. the true signal. Hence in applications, the decrease in energy of the residual 

is often treated as a measure of performance improvement. Proposition [2] gives a sufficient 


April 8, 2015 


DRAFT 

















13 


condition for decrease in the energy of the residual matrix obtained by MMV-FACS over the z lh 
participating algorithm. 


Proposition 2. For a K-sparse matrix X, let R = B — A r At.B and R, = B — Af, At.B repre- 

\/l + £ 1+ 


sent the residue matrix of MMV-FACS and i th Algorithm respectively. Assume that 


JR+K 


8r+k + 3£) < 


1-2 5 


R+K 


fli \/1 — $R+K 

Proof: We have, 


— £ is satisfied then we have, ||R|L < ||R, 


i\\F- 


IRI 


Using 


B - A r At B 


AX + W - A^At (AX + W) (y B = AX + W) 

T F 


A ,X, . 


A rc X tCj : - A t A^A tc X tCj : + W - A t A\AV 


A tc X tC): + W - A t At, (A t X t: + A tc X ts: + W 

" f (-' a W=i) 


< 

< 


A j-cXj-c . Ay-A^ A«y- C X^ c . 


W - A t At,'W 


A| C X| C . 


A r\f Xr\r,: 


|W|| F (Using Lemma [2]) 


IWj 


< v'T+w 

we have, 


X 


T c ,: 


+ II w| 


r = supp(x)) 


v|T\t|<x & s K <s R+K 


|R|| f < y/T+W ( ? + Sr+K ||Xr«= 


1-5 


R+K 


If+ 1-5 


IWI 


R+K 


+ II w| 


< 


< 


\/l + 8 


R+K 


1 — 8r+k 

y/l + 8 r+k 


1-5 


R+K 


1 + 8 r+ k + 2£ H— R+ =£ ) ||X r c : || F 

y/ 1 + 5 R+k 

(1 + 8 r+ k + 3£) ||Xrc >: || F . 


(31) 
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|Ri|| F = 


B — A~-At B 

% Ti 


AX + W — At (AX + W) (y B = AX + W) 


A^X^ . 

I i lit' 


Ti 

AfcXfc- + W - A^At ('A, X /; : + A+ w 


A^X^ . + W - A r At A^X^ . - A, r At 'W 

' i 'i > • '* Ti 'i 'i > • !% Ti 


> 


A tc X* . - A-fAt A FC X FC . 

'i 'i »• '* 7i 'i 'i ’• 


W - A.+-At W 

7i 7i 


(Using reverse triangle inequality) 


A T\T^T\ti,-. A ^ A ^ A r\7i X r\7i,: 


> 1 - 


J-R+AT 


A T\T X TV7Y 


> 


1 — Sr+K 
(Using Lemma [2] & 5 2 k < <5r+a 
1 — 2<5 F+F 


F 

-||W|| F 


W- At-A t'W 


1 — ^R+K 

1 — ^R+K 




R+K 


X 


TV 


IWI 


{■:\T\T\<K & 5 k <8 


R+K) 


C l|Xp=, : 


If • 


J>li \/1 — <tfM-A 

From (OTb and (l32l) we get a sufficient condition for ||R|| F < ||R, || F as 

\/l + 8 R+K 


(32) 


■(1 + <5r+a' + 3£) < 


1 — 25 f+ a' 
r h \T^~~5 r+k 


-c 


(33) 


1 — 8 R+K 

Thus, if (l33l) is satisfied, MMV-FACS produces a smaller residual matrix (in the Frobenius norm 
sense) than that of the z th participating algorithm. ■ 


C. Average Case Analysis 

Intuitively, we expect multiple measurement vector problem to perform better than the single 
measurement vector case. However, if each measurement vector is the same, i.e., in the worst 
case, we have = c, V i = 1,..., L, then we do not have any additional information on 
X than that provided by a single measurement vector So far we have carried out only 
the worst case analysis, i.e., conditions under which the algorithm is able to recover any joint 
sparse matrix X. This approach does not provide insight into the superiority of sparse signal 
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reconstruction with multiple measurement vectors compared to the single measurement vector 
case. 

To notice a performance gain with multiple measurement vectors, next we proceed with an 
average case analysis. Here we impose a probability model on the K sparse X as suggested by 
Remi et al. ifTTIl . In particular, on the support-set T, we impose that X 7 -. = £$, where £ is a 
K x K diagonal matrix with positive diagonal entries and $ is a I\ x L random matrix with 
i.i.d. Gaussian entries. Our goal is to show that, under this signal model, the typical behaviour 
of MMV-FACS is better than in the worst case. 


Theorem 2. Consider the MMV-FACS setup discussed in Section [ 777] Assume a Gaussian signal 


model, i.e., X 7 -. : = £$, where £ is a K x K diagonal matrix with positive diagonal entries and 


$ is a K x L random matrix with i.i.d. Gaussian entries. Let e, denote a ITI x 1 vector with a 


‘1’ in the i th coordinate and ‘0’ elsewhere. Let q 
and 


min 

ie(Tnr) 


e fApW 


2 


+ max 
je(T\T) 


ej Ap W 


2 


7 


min 

ie(Tnr) 


ef Ap A r £ 


2 


max 

J6(T\T) 


ej Ap A r £ 


2 


V 

c 2 (L) 


min 

ej Ap A r £ 

+ max 

ej Ap A r £ 


iG(rnr) 


2 je(r\T) 


2 


where C 2 (L) = E \\Z || 2 with Z = (Zi,..., Zj) being a vector of independent standard normal 


variables. Assume that min 

ie(Tnr) 


ef Ap A r £ 


2 


max 

je(r\T) 


ejApA r £ 


> . Let 0 denote 

Ci\L) 


the event that MMV-FACS picks all correct indices from the union-set T. Then, we have, 


P(0) > 1-/Texp(-2A 2 (L) 7 2 ), 


/ T(C±±) 
where A 2 (L) — I .. 2 — 

V r (f) 

Proof: We have, 


L ... . 

—, T(-) denotes the Gamma function. 


P(Q) — P I min 
' ie(rnr) 


= P I min 

^je(Tnr) 

> P I min 
\ie(Tnr) 

> max 

je(r\r) 


efA*B 


> max 
2 ie(r\r) 


ejApB 


efA^(A r X T , : +W) 
e f ApA r X ri . 


> max 
2 je(r\T) 


e I A|, W 


e[A{,(A r X r , : +W) 


e J Ap A r X ri . 


+ 


ejApW 
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(Using reverse triangle inequality and triangle inequality respectively) 


= P ( min 

Vie(rnr) 

> min 

jg(Tnr) 

= P ( min 

jg(Tnr) 


e f Ap A7- X T) . 


— max 
2 je(r\r) 


ej Ap A r X Ti: 


e fApW 


+ max 
2 je(r\r) 


ejApW 


e f ApA r X Ti . 


— max 
2 je(r\T) 


e jApA r X T) : 


> rj 


= l-p 

min 

\ie(Tnr) 

ef A[ A t X T: . 

>1 — p 

min 

\ie(T nr) 

e f ApA r X Tj . 

- p ( 

max 

yi6(r\T) 

e J Ap A r X T): 


— max 
2 je(r\r) 


< C 


e jApA r X ri . 


< r] 


>C-rj). 


Now, let us derive an upper bound for P min 

V* e ( rnr ) 

concentration of measure results in lUTTl . we set 

C = (1 - €1 )C 2 (L) min 
*e(Tnr) 


ef AtA r X r : 


e I ApA- 7 - E 


(34) 

< C ). Influenced by the 

(35) 


where 0 < ei < 1. 

Using (5.5) in 0, we get, 


P ( min 

\ie(mr) 


Va*a t x t , 


2 



< \T\exp(-A 2 (L)e 2 1 ). 


(36) 


To bound the second probability, consider 


P I max 
vje(ryr) 


e^ApA r S<I) 



P max 

e^ApArE® 

1 je(r\r) 



>(C-v) 


Co(L) max 
je(r\T) 


Co(L) max 

je(r\T) 


ejApA r E 

ejApA r E 


\ 


2/ 


Let 


1 + e 2 


C-T) 


Co ( L) max 
je(r \T) 


ej Ap A r E 


2 


(37) 


Using equation (5.3) in [fTTll 


P | max 

pe(r\T) 


ejApA r E$ 


> (1 + e 2 )C 2 (L) max 
2 J 6 (T\T) 


ejAfA r E 


< 17] exp (~A 2 (L)e 2 ) . 


(38) 
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For the above inequality to hold, it is required that e 2 > 0. By setting e 2 = ei, and using ([331) 
and (137b . we get 


(1 - 6i)C 2 (L)min ie(rnr) 

efA*A 

r s 

-T] 

2 

(7 2 (L) max 
jg(t \T) 

ej Ap A r E 

2 


Now, solving for e, we get 


ei = 


min 

ie(Tnr) 


Ap A r E 


2 


max 
i e(r\r) 


ejApA r E 


2 


c 2 {l) 


min 

jg(Tnr) 


e f Ap A r E 


+ max 
2 je(r\r) 


ejApA r E 


2 


Clearly ei < 1 and by the assumption in the theorem e\ > 0. Hence we have 0 < e\ < 1. Also, 


note that 7 = e 1 . Substituting (l36b and (l38b in (l34b . we get 


P(0) > 1 - K exp(—2A 2 (L) 7 2 ). 


Since A 2 (L) « —, the probability that MMV-FACS selects all correct indices from the union set 
increases as L increases. Thus, more than one measurement vector improves the performance. 

V. Numerical Experiments and Results 

We conducted numerical experiments using synthetic data and real signals to evaluate the 
performance of MMV-FACS. The performance is evaluated using ASRER which is defined as 

E ^trials 

.7=1 


ASRER = 


ry u 

<7=1 




2 

\ ''^'trials 
2 ^j=i 

x, - x, 

F 


(39) 


where X ? and Xj denote the actual and reconstructed jointly sparse signal matrix in the j th trial 
respectively, and n tri , a i s denotes the total number of trials. 


A. Synthetic Sparse Signals 


For noisy measurement simulations, we define the SMNR as 

| 2 i 

1 2 J 


a e{||x«||;} 


SMNR = 9 , 

E{||wW|| 2 } 

where E{-} denotes the mathematical expectation operator. The simulation set-up is described 
below. 
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1) Experimental Setup: Following steps are involved in the simulation: 

1) Generate elements of A MxN independently from jj) and normali z e each column 
norm to unity. 

ii) Choose K non-zero locations uniformly at random from the set {1, 2,, N} and fill those 
K rows of X based on the choice of signal characteristics: 

a) Gaussian sparse signal matrix: Non-zero values independently from J\f( 0, 1). 

b) Rademacher sparse signal matrix: Non-zero values are set to +1 or -1 with probability 

i 

2 ' 

Remaining N — K rows of X are set to zero. 

iii) The MMV measurement matrix B is computed as B = AX+W, where the columns of W, 
wW’s are independent and their elements are i.i.d. as Gaussian with variance determined 
from the specified SMNR. 

iv) Apply the MMV sparse recovery method. 

v) Repeat steps i-iv, S times. 

vi) Find ASRER using (l39l) . 

2) Results and Discussions: We used M-OMP, M-SP, M-BPDN lfl8l . and M-FOCUSS CD as 
the participating algorithms in MMV-FACS. The software code for M-BPDN was taken from 
SPGL1 software package |fl9l . Since M-FOCUSS and M-BPDN algorithms may not yield an 
exact /v'-sparsc solution, we estimate the support-set as the indices of the K rows with largest £2 
norm. We fixed the sparse signal dimension N = 500 and sparsity level K = 20 in the simulation 
the result were calculated by averaging over 1, 000 trials (S = 1,000). We use an oracle estimator 
for performance benchmarking. The oracle estimator is aware of the true support-set and finds 
the non-zero entries of the sparse matrix by solving LS. 

The empirical performance of MMV reconstruction algorithms for different values of M is 
shown in Fig. HI The simulation parameters are L = 20, SMNR= 20 dB and X is chosen as 
Gaussian sparse signal matrix. For M = 35, MMV-FACS (M-BPDN,M-FOCUSS) gave 10.67 
dB and 4.27 dB improvement over M-BPDN and M-FOCUSS respectively. 

Fig. [2] depicts the performance of Rademacher sparse signal matrix for different values of M 
where we set L = 20 and SMNR= 20 dB. We again observe similar performance improvement 
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(b) 



Number of Measurements (M) 


Number of Measurements (M) 


(c) 


(d) 


Fig. 1. Performance of MMV-FACS, averaged over 1,000 trials, for Gaussian sparse signal matrices with SMNR = 20 dB. 
Sparse signal dimension N = 500, sparsity level K = 20, and number of measurement vectors L = 20. 


as in the case of Gaussian sparse signal matrix. For example, for M = 35, MMV-FACS(M- 
OMP,M-BPDN) showed 7.56 dB and 4.32 dB over M-OMP and M-BPDN respectively. 

A comparison of MMV reconstruction techniques is shown in Fig.[3jbr Gaussian sparse signal 
matrix for different values of L where we set M — 50 and SMNR= 20 dB. It may be observed 
that MMV-FACS gave a significant performance improvement over the participating algorithms. 
Specifically, MMV-FACS(M-OMP,M-SP) improved the performance by 5.77 dB and 4.94 dB 
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(a) 


(b) 



Number of Measurements (M) 


Number of Measurements (M) 


(c) 


(d) 


Fig. 2. Performance of MMV-FACS, averaged over 1000 trials, for Rademacher sparse signal matrices with SMNR = 20 dB. 
Sparse signal dimension N = 500, sparsity level K = 20 and number of measurement vectors L = 20. 


over M-OMP and M-SP respectively. 

To show the dependency of recovery performance on SMNR, we conducted simulations for 
different values of SMNR. Fig. 0] illustrates the performance for Gaussian sparse signal matrix 
where L — 10 and M = 45. An additional ASRER improvement of 2.51 dB and 2.08 dB were 
achieved as compared to M-OMP and M-FOCUSS respectively for SMNR= 10 dB. This shows 
the robustness of MMV-FACS to noisy measurements. 
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Number of Measurement Vectors (L) Number of Measurement Vectors (L) 


(a) 


(b) 


Fig. 3. Performance of MMV-FACS, averaged over 1000 trials, for Gaussian sparse signal matrices with SMNR = 20 dB. 
Sparse signal dimension N = 500, sparsity level K = 20, and number of measurements M = 50. 


From the above simulation results it can be seen that MMV-FACS improved the sparse signal 
recovery compared to participating algorithms. 

3) Reproducible Research: We provide necessary Matlab codes to reproduce all the figures, 
publicly downloadable from http://www.ece.iisc.ernet.in/~ssplab/Public/MMVFACS.tar.gz. 

B. Real Compressible Signals 

To evaluate the performance of MMV-FACS on compressible signals and real world data, we 
used the data set ‘ 05091 .daf from MIT-BIH Atrial Fibrillation Database [f20ll . The recording 
is of 10 hours in duration, and contains two ECG signals each sampled at 250 samples per 
second with 12-bit resolution over a range of ±10 millivolts. We selected the first 250 time 
points of the recording as the data set used in our experiment. We used a randomly generated 
Gaussian sensing matrix of size M x 250, with different values of M in the experiment. We 
assumed sparsity level K = 50 and used M-OMP and M-SP as the participating algorithms. The 
reconstruction results are shown in Fig. [5] 

Similar to synthetic signals, MMV-FACS shows a better ASRER compared to the participating 
algorithms M-OMP and M-SP. This demonstrates the advantage of MMV-FACS in real-life 
applications, requiring fewer measurement samples to yield an approximate reconstruction. 
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(c) (d) 


Fig. 4. Performance of MMV-FACS, averaged over 1000 trials, for Gaussian sparse signal matrices with SMNR = 20 dB. 
Sparse signal dimension N = 500, sparsity level K = 20, and number of measurements M = 45, and number of measurement 
vectors L = 10. 


VI. Conclusions 

In this paper, we extended FACS to the MMV case and showed that MMV-FACS improves 
sparse signal matrix reconstruction. Using RIP, we theoretically analysed the proposed scheme 
and derived sufficient conditions for the performance improvement over the participating algo¬ 
rithm. Using Monte-Carlo simulations, we showed the performance improvement of the proposed 
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Fig. 5. Real Compressible signals: Performance of MMV-FACS for 2-channel ECG signals from MIT-BIH Atrial Fibrillation 
Database (201. 


scheme over the participating methods. Though this paper discusses only the extension of FACS 
for MMV problem, a similar approach can be used to extend the other fusion algorithms 
developed by Ambat et al. l!T5l . OTt ll22ll . 


Appendix 

Proof of Lemma □ 

The proof is inspired by Proposition 3.5 by Needell and Tropp lfl6l . 

Define set S as the convex combination of all matrices which are R + K sparse and have unit 
Frobenius norm. 

S = convI^X : ||X|| 0 < R + K, ||X|| F = lj 
Using the relation ||AX|| f < + 5r+k ||X|| f , we get. 


Define 


argmax || AX IIf <yiT &R+K ■ 
xgs 


Q 



IXIL + 


VR + K 


|X|| 2il < 
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The lemma essentially claims that 


argmax || AX|| F < argmax ||AX|| F . 

xgQ xes 

To prove this, it is sufficient to ensure that Q C S. 

Consider a matrix X e Q. Partition the support of X into sets of size R + K. Let set X 0 
contain the indices of the R + K rows of X which have largest row f? 2 -norm, breaking ties 
lexicographically. Let set X\ contain the indices of the next largest (row f^-norm) R + K rows 
and so on. The final block Ij may have lesser than R + K components. This partition gives rise 
to the following decomposition: 

j J 

X = X |/ 0 + ^ x|j. = AoYo + ^ XjYj, 

3 =1 3 =1 

where Xj = ||X|/J F and Yj = A~' X|/.. 

By construction each matrix Y j belongs to S because it is R+ K sparse and has unit Frobenius 
norm. We will show that JY \. < i. This implies that X can be written as a convex combination 
of matrices from the set S. As a result X G S. Therefore, Q C S. 

Fix some j in the range {1,2,..., J}. Then, Ij contains at most R + K elements and 2}| 
contains exactly R + K elements. Therefore, 


Aj = Wi.\\ F <y/R+K\\X\i. 


, lln <VR + K - 

J 2,oo V R + R 


XI 


j- 1 112,1 ' 


The last inequality holds because the row £ 2 -norm of X on the set / ; _ | dominates its largest 
row £ 2 -norm in Ij. Summing these relations, we get 

J 1 n 1 




, < 


XIAiv.lly 


< 


VR + K 


X 


\/R + K “ 

3 =1 v 3 = 1 

Also, we have A 0 = ||X|/ 0 || F < ||X|| F . Since X G Q, we conclude that 

j 


2,1 • 


X a j^ 


3=0 


|X|L H — T =^= ||X| 

VrTk 


2,1 


< 1 . 


Proof of Lemma [2] 

Let yW denote the / th column of matrix Y and r'Y denote the z lh column of matrix R, 
i = 1, 2,... L. Then we have from Lemma 2 of 


s \Ti\+m\ \ 
1 - Vil+ITal/ 




< 

2 — 


.(*) I 


< ||y 


(i)| 
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Summing the above relation, we obtain 


fyn\±M\ A 
- Vil+ir 2 | J “ 


r(0 


L 

I un 

i= 1 


.(<) 


L 

\UT. 

%— 1 


r(0 


2 

2 ' 


Equivalently, we have, 



Vil+ir 2 | \ 

1 - ViI+|T 2 |/ 


|y|| f < 


|R|| f < ||Y|| f . 
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